package scu.maqiang.mesh;

import java.util.HashSet;

import scu.maqiang.numeric.CG;
import scu.maqiang.numeric.MVO;
import scu.maqiang.numeric.ParamCheck;
import scu.maqiang.numeric.UnorderIntArray;

import static java.lang.Math.*;
import static java.lang.Math.sin;

/**
 * 
 * 二维四结点四边形网格
 * 
 * @author 马强
 *
 */
public class Mesh2Q4 extends Mesh {

	public Mesh2Q4() {
		super();
		nDim = 2;
		nPerElement = 4;
		nPerBoundary = 2;
		nEdgePerElement = 4;
		nEdgePerBoundary = 1;
		tecplotType = "FEQUADRILATERAL";
		nBoundaryPerElement = 4;
		meditMeshType = new String[]{"Vertices", "Quadrilaterals", "Edges"};
		gmshType = new String[]{"ST", "VT", "TT"};
		gmshElementType = 3;
	}

	public Mesh2Q4 square2D(int xNum, int yNum) {
		int nv = (xNum + 1) * (yNum + 1);
		int nt = xNum * yNum;
		int nb = 2 * (xNum + yNum);
		initMesh0(nv, nt, nb);
		int count = 0;
		for (int i = 0; i < yNum + 1; i++) {
			for (int j = 0; j < xNum + 1; j++) {
				nodes[count++] = new double[]{(j + 0.0) / xNum, (i + 0.0) / yNum};
			}
		}
		createSquareElements(xNum, yNum);
		createSquareBoundarys(xNum, yNum);
		return this;
	}

	public Mesh2Q4 quadrilateral(double[] QuadX, double[] QuadY, int xNum, int yNum) {

		int n, i, j, xNodes, yNodes, ti, tj, tk, tl;
		double xi, et, x, y;
		double[] h = new double[4];
		double[] dXi = new double[4];
		double[] dEt = new double[4];
		double[] hh = new double[8];
		double[] ddXi = new double[8];
		double[] ddEt = new double[8];
		n = QuadX.length;
		xNodes = xNum + 1;
		yNodes = yNum + 1;
		int nv = xNodes * yNodes;
		int nt = xNum * yNum;
		int nb = 2 * (xNum + yNum);
		initMesh0(nv, nt, nb);
		int count = 0;
		if (n == 4) {
			count = 0;
			for (j = 0; j < yNodes; j++) {
				for (i = 0; i < xNodes; i++) {
					xi = 2 * (i + 0.0) / xNum - 1.0;
					et = 2 * (j + 0.0) / yNum - 1.0;
					CG.Q4ShapeFunction(xi, et, h, dXi, dEt);
					x = MVO.dot_product(h, QuadX);
					y = MVO.dot_product(h, QuadY);
					nodes[count++] = new double[] {x, y};
				}
			}
		} else if (n == 8) {
			count = 0;
			for (j = 0; j < yNodes; j++) {
				for (i = 0; i < xNodes; i++) {
					xi = 2 * (i + 0.0) / xNum - 1.0;
					et = 2 * (j + 0.0) / yNum - 1.0;
					CG.Q8ShapeFunction(xi, et, hh, ddXi, ddEt);
					x = MVO.dot_product(hh, QuadX);
					y = MVO.dot_product(hh, QuadY);
					nodes[count++] = new double[] {x, y};
				}
			}
		} else
			throw new IllegalArgumentException("The length of the array must be 4 or 8!");

		count = 0;
		for (j = 0; j < yNum; j++) {
			for (i = 0; i < xNum; i++) {
				ti = i + j * xNodes;
				tj = ti + 1;
				tk = tj + xNodes;
				tl = tk - 1;
				elements[count++] = new int[] {ti, tj, tk, tl};
			}
		}

		createSquareBoundarys(xNum, yNum);
		return this;
	}


	public Mesh2Q4 ring(double ra, double rb, int nR, int nL) {
		double[] theta = MVO.linspace(0, 2 * PI, nR + 1);

		int nv = (nL + 1) * nR;
		int nt = nR * nL;
		int nb = 2 * nR;
		initMesh0(nv, nt, nb);
		double dr = (rb - ra) / nL;
		double x, y;
		int count = 0;
		for (int i = 0; i < nL + 1; i++) {
			for (int j = 0; j < nR; j++) {
				x = (ra + i * dr) * cos(theta[j]);
				y = (ra + i * dr) * sin(theta[j]);
				// System.out.println(MVector.toString(x, y));
				nodes[i] = new double[] {x, y};
			}
		}

		count = 0;
		int L1, L2, L3, L4;
		for (int i = 0; i < nL; i++) {
			for (int j = 0; j < nR - 1; j++) {
				L1 = j + i * nR;
				L2 = L1 + nR;
				L3 = L2 + 1;
				L4 = L1 + 1;
				elements[count++] = new int[] {L1, L2, L3, L4};
			}
			L1 = i * nR + nR - 1;
			L2 = L1 + nR;
			L3 = (i + 1) * nR;
			L4 = i * nR;
			elements[count++] = new int[] {L1, L2, L3, L4};
		}

		count = 0;
		for (int i = 0; i < nR - 1; i++) {
			L1 = i;
			L2 = i + 1;
			boundaries[count] = new int[] {i, i + 1};
			boundaryLabel[count++] = 1;
		}
		L1 = nR - 1;
		L2 = 0;
		boundaries[count] = new int[] {L1, L2};
		boundaryLabel[count++]= 1;

		int local = nL * nR;
		for (int i = 0; i < nR - 1; i++) {
			L1 = i + nL * nR;
			L2 = L1 + 1;
			boundaries[count] = new int[] {local, local + 1};
			boundaryLabel[count++] = 2;
		}

		L1 = nL * nR + nR - 1;
		L2 = nL * nR;
		boundaries[count] = new int[] {L1, L2};
		boundaryLabel[count] = 2;
		return this;
	}

	public Mesh2Q4 fromRectangleNodes(double[][] nodes, int xNum, int yNum) {
		nv = nodes.length;
		ParamCheck.checkEqual(nv, (xNum + 1) * (yNum + 1));
		nt = xNum * yNum;
		nb = 2 * (xNum + yNum);
		initMesh0(nv, nt, nb);
		for(int i = 0; i < nv; i++) {
			nodes[i] = nodes[i].clone();
		}
		createSquareElements(xNum, yNum);
		createSquareBoundarys(xNum, yNum);
		return this;
	}

	public Mesh2Q4 fan(double ra, double rb, double rc, double rd, double theta1, double theta2,
					   int rNum, int thetaNum, T3Type type) {
		nv = (rNum + 1) * (thetaNum + 1);
		nt = rNum * thetaNum;
		nb = 2 * (rNum + thetaNum);
		initMesh0(nv, nt, nb);
		double rab = rb - ra;
		double rdc = rc - rd;
		double dr0 = rab / rNum;
		double dr1 = rdc / rNum;

		int count = 0;
		for (int i = 0; i < thetaNum + 1; i++) {
			double dr = dr0 + i * (dr1 - dr0) / thetaNum;
			double R0 = ra + i * (rd - ra) / thetaNum;
			double theta = theta1 + i * (theta2 - theta1) / thetaNum;
			for (int j = 0; j < rNum + 1; j++) {
				double R = R0 + dr * j;
				nodes[count++] = new double[] { R * cos(theta), R * sin(theta)};
			}
		}

		createSquareElements(rNum, thetaNum);
		createSquareBoundarys(rNum, thetaNum);
		return this;
	}

	private void createSquareElements(int xNum, int yNum) {
		int ii, jj, kk, ll;
		int count = 0;
		for (int i = 0; i < yNum; i++) {
			for (int j = 0; j < xNum; j++) {
				ii = i * (xNum + 1) + j;
				jj = ii + 1;
				kk = jj + xNum + 1;
				ll = ii + xNum + 1;
				elements[count++] = new int[] {ii, jj, kk, ll};
			}
		}
	}

	private void createSquareBoundarys(int xNum, int yNum) {
		int count = 0;
		for (int i = 0; i < xNum; i++) {
			boundaries[count] = new int[] {i, i + 1};
			boundaryLabel[count++] = 1;
		}
		int local = xNum;
		for (int i = 0; i < yNum; i++) {
			boundaries[count] = new int[] {local, local + xNum + 1};
			boundaryLabel[count++] = 2;
			local += xNum + 1;
		}

		local = xNum + yNum * (xNum + 1);
		for (int i = 0; i < xNum; i++) {
			boundaries[count] = new int[] {local, local - 1};
			boundaryLabel[count++] = 3;
			local--;
		}

		local = yNum * (xNum + 1);
		for (int i = 0; i < yNum; i++) {
			boundaries[count] = new int[] {local, local - (xNum + 1)};
			boundaryLabel[count++] = 4;
			local -= xNum + 1;
		}
	}

//	public Mesh2DQ4 add(Mesh2DQ4 meshB) {
//	    int[] newNodesB, newBdB, ABd;
//	    double tol = 1.0e-3, xA, yA, xB, yB;
//	    int i, j, nvA, nvB, nConnected, n, nBd;
//	    int[] idx = new int[2];
//	    int[] idx2 = new int[2];
//	    int[] temp = new int[2];
//	    newNodesB = new int[meshB.nv];
//	    newBdB = new int[meshB.nb];
//	    ABd = new int[nb];
//	    Arrays.fill(newNodesB, -1);
//	    Arrays.fill(newBdB, -1);
//	    Arrays.fill(ABd, -1);
//	    nvA = nv;
//	    nvB = meshB.nv;
//	    for(i = 0; i < nvA; i++) {
//	        xA = nodes.get(i).coord[0];
//	        yA = nodes.get(i).coord[1];
//	        for(j = 0; j < nvB; j++) {
//	            xB = meshB.nodes.get(j).coord[0];
//	            yB = meshB.nodes.get(j).coord[1];
//	            if (abs(xA - xB) < tol && abs(yA - yB) < tol) {
//	                newNodesB[j] = i;
//	                break;
//	            }
//	        }
//	    }
//
//	    for(i = 0; i < meshB.nb; i++) {
//	        idx = meshB.boundarys.get(i).idx;
//	        temp[0] = newNodesB[idx[0]];
//	        temp[1] = newNodesB[idx[1]];
//	        if (temp[0] > -1 && temp[1] > -1) {
//	            newBdB[i] = 1;
//	            for(j = 1; j < nb; j++) {
//	                idx2 = boundarys.get(j).idx;
//	                if ((idx2[0] == temp[0] && idx2[1] == temp[1]) ||
//	                    (idx2[1] == temp[0] && idx2[1] == temp[0])) {
//	                       ABd[j] = 1;
//	                       break;
//	                }
//	            }
//	        }
//	    }
//
//	    //nConnected = (int)(IntStream.of(newNodesB).filter(iidx -> iidx != -1).count());
//	    //nBd = (int)(IntStream.of(newBdB).filter(iidx -> iidx != -1).count());
//	    nConnected = 0;
//	    n = nvA;
//	    for(i = 0; i < nvB; i++) {
//	        //���B�еĽ��û����A�еĽ����ԣ�����������B�еĽ�����µ������еĽ����
//	        if (newNodesB[i] == -1) {
//	            newNodesB[i] = n;
//	            n++;
//	        }else {
//	        	//����Ѿ���ԣ����������ԵĽ��ĸ���
//	        	nConnected++;
//	        }
//	            
//	    }
//	    
//	    nBd = 0;
//	    n = nb;
//	    for(i = 0; i < meshB.nb; i++) {
//	    	if (newBdB[i] == -1) {
//	            newBdB[i] = n;
//	            n++;
//	    	}
//	        else {
//	            nBd++;
//	        }
//	    }
//	    
//	    nv = nvA + nvB - nConnected;
//	    nt = nt + meshB.nt;
//	    nb = nb + meshB.nb - nBd;
//	    
//	    for(i = 0; i < nvB; i++) {
//	    	if (newNodesB[i] > nvA - 1) {
//	    		Node node = meshB.nodes.get(i);
//	    		nodes.add(new Node3D((Node3D)node));
//	    	}
//	    }
//	    
//	    for(i = 0; i < meshB.nt; i++) {
//	    	Quadrilateral4 ele = (Quadrilateral4)meshB.elements.get(i);
//	    	int[] newEle = new int[4];
//	    	for(int kk = 0; kk < 4; kk++) {
//	    		newEle[kk] = newNodesB[ele.idx[kk]];
//	    	}
//	    	elements.add(new Quadrilateral4(newEle, ele.label));	    	
//	    }
//	    
//	    for(i = 0; i < boundarys.size(); i++) {
//	    	//System.out.println(i + "\t" + nb + "\t" + boundarys.size());
//	    	if(ABd[i] == 1) {
//	    		boundarys.get(i).label = 0;
//	    	}
//	    }
//	    
//	    for(i = 0; i < meshB.nb; i++) {
//	    	if (newBdB[i] != 1) {
//	    		Element li = meshB.boundarys.get(i);
//	    		int[] newEle = new int[2];
//	    		for(int kk = 0; kk < 2; kk++) {
//	    			newEle[kk] = newNodesB[li.idx[kk]];
//	    		}
//	    		boundarys.add(new Line2(newEle, li.label));
//	    	}
//	    }
//		return this;
//	}
//	
//	public Mesh2DQ4 add(Mesh2DQ4... mesh) {
//		for(int i = 0; i < mesh.length; i++) {
//			this.add(mesh[i]);
//		}
//		return this;
//	}

	@Override
	public double[] getinteriorPointInElement(double[][] coord) {
		double[] x = new double[4];
		double[] y = new double[4];
		for (int i = 0; i < 4; i++) {
			x[i] = coord[i][0];
			y[i] = coord[i][1];
		}
		return new double[] { MVO.average(x), MVO.average(y) };
	}

	@Override
	public double[] getinteriorPointOnBoundary(double[][] coord) {
		double[] x = new double[2];
		double[] y = new double[2];
		for (int i = 0; i < 2; i++) {
			x[i] = coord[i][0];
			y[i] = coord[i][1];
		}
		return new double[] { 0.5 * MVO.sum(x), 0.5 * MVO.sum(y) };
	}

	@Override
	public boolean checkMesh(double[][] coord) {
		return CG.area2DQ4(coord) > 0.0;
	}

	@Override
	public int[][] getBoundarysFromElement(int i) {
		int[][] boundarys = new int[nBoundaryPerElement][nPerBoundary];
		int[] ele = elements[i];
		boundarys[0][0] = ele[0];
		boundarys[0][1] = ele[1];
		boundarys[1][0] = ele[1];
		boundarys[1][1] = ele[2];
		boundarys[2][0] = ele[2];
		boundarys[2][1] = ele[3];
		boundarys[3][0] = ele[3];
		boundarys[3][1] = ele[0];
		return boundarys;
	}

	public int[] signsInEdges(int[] ele) {
		int[] signs = new int[4];
		signs[0] = ele[0] < ele[1] ? 1:-1;
		signs[1] = ele[1] < ele[2] ? 1:-1;
		signs[2] = ele[2] < ele[3] ? 1:-1;
		signs[3] = ele[3] < ele[0] ? 1:-1;
		return signs;
	}

	@Override
	public void setBoundaries(HashSet<UnorderIntArray> set) {
		this.nb = set.size();
		boundaries = new int[this.nb][];
		boundaryLabel = new int[this.nb];
		int count = 0;
		for (UnorderIntArray ele : set) {
			boundaries[count++] = ele.getIntSet();
		}
	}

	@Override
	public int getTecplotNt() {
		return nt;
	}

	@Override
	public int[][] getEdgesFromElement(int i) {
		int[] ele = elements[i];
		int[][] result = { { ele[0], ele[1] }, { ele[1], ele[2] }, { ele[2], ele[3] }, { ele[3], ele[0] } };
		return result;
	}

}
